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K*" Abstract 

o 

We consider a finite element method which couples the continuous Galerkin method 
away from internal and boundary layers with a discontinuous Galerkin method in the 
vicinity of layers. We prove that this consistent method is stable in the streamline 
I— I diffusion norm if the convection field flows non-characteristically from the region of 

the continuous Galerkin to the region of the discontinuous Galerkin method. The 
stability properties of the coupled method are illustrated by numerical experiments. 

1 Introduction 

B 

I— I It is well known that the standard continuous Galerkin (cG) finite element method exhibits 

I poor stability properties for singularly perturbed problems. In presence of sharp boundary 

^ or interior layers, non physical oscillations pollute the numerical approximation throughout 

the solution domain, and thus stabilisation techniques need to be employed; see [13] for 
a survey. The discontinuous Galerkin (dG) method offers a framework for the design of 
finite element methods with good stability properties; see, e.g., |8] for a survey of their 
7—1 development. This is achieved by relaxing the continuity requirements at the inter-element 

boundaries, where appropriate upwinded numerical fluxes can be employed. One such dG 
method is the Interior Penalty (IP) method considered herein. 

However dG methods require more degrees of freedom (DOFs) compared to the stan- 
dard cG method. In this work, we investigate if some of these additional DOFs can be 
^ removed without affecting the stability properties of the dG method. The motivation of 

this work is twofold: firstly we wish to improve our understanding of the well-established 
dG methods in response to the criticism in the increased number of DOFs; secondly we 
seek insight into the design of more efficient finite element methods, allowing to combine 
advantages the dG framework with smaller approximation spaces. 

Conceptually we envisage in this work a finite element space that lies between the 
standard continuous and discontinuous Galerkin spaces. We construct such a space by 
applying standard continuous elements away from any boundary or internal layers (called 
the cG region) and discontinuous elements in the region of such layers (called the dG 
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region). Therefore we call this method the continuous discontinuous Galerkin (cdG) finite 
element method. 

In the broader sense, methods of this type have been proposed before by Becker, 
Burman, Hansbo and Larson [Ij, Perugia and Schotzau [13], and Dawson and Proft [10]. 
In [3] a globally reduced discontinuous Galerkin method is studied, whose approximation 
space consists of the continuous piecewise linear functions enriched with piecewise constant 
functions. In |13j . |10| the local discontinuous Galerkin method [9] has been used on the 
discontinuous region with transmission conditions on the subset of interelement boundaries 
where continuous and discontinuous elements meet. Our approach is different to |13j . 
[TU] as we impose no such conditions beyond those already imposed by the dG method. 
Removing the transmission conditions makes the approach more natural, but introduces 
some difficulties as possible over- and undershoots at the interface of the cG and dG regions 
must be controlled. Indeed, the analysis is limited to the case where the cdG interface 
is non-characteristic so that jumps across the interface are controlled by the convection 
term. This allows us to derive a rigorous stabilty bound for the cdG approximation. On 
the cG region we split the numerical solution into an approximation to the hyperbolic 
problem which is weakly dependent on the diffusion coefficient e and an approximation to 
the remaining part, which depends more strongly on e but is small in size. Stability of the 
approximation on the dG region is shown using the approach of Buffa, Hughes and Sangali 
[5] with some extensions in the manner of Ayuso and Marini [2]. To our knowledge this 
paper is the first to present a stability result for the proposed cdG method. 

Comparisons of our cdG method with various dG and cG methods have been under- 
taken by Cangiani, Georgoulis and Jensen ffj and Devloo, Forti and Gomes |llj ; see also 
[6] which studies the cdG approximation as the limit of a dG approximation as the jump 
penalization on interelement boundaries tends to infinity. 

The remainder of this work is organized as follows. In Section [2] we formulate the prob- 
lem, introduce notation, and explain the assumptions we make. Control on the continuous 
region is considered in several stages in Section [3] and Section |4| Section [5] presents an 
inf-sup bound on the discontinuous region. The main result of the paper showing control 
independent of the perturbation parameter e and mesh size h follows in Section [6} The 
theoretical results are illustrated in Section [7] by a numerical experiment. 

2 The Model Problem and Notation 

Let i7 be a bounded Lipschitz domain in M*^. We introduce the model advection-diffusion- 
reaction (ADR) problem 



with constant diffusion coefficient < e < Smax, b G [W^'°°{Q)f, c G L°°{Q) and / G 
L'^(Q). For < e ^ 1 the solution to this problem typically exhibits boundary or interior 
layers. 

Unless otherwise stated, we define C throughout as a positive constant, independent 
of £ and the finite element approximation space and which may be redefined from line to 
line. C may depend on Emax- By a < b we mean a < Cb. 

We consider Th to be subdivision of into non-overlapping shape regular simplices or 
hyper-cubes E, which we shall refer to as the triangulation. Denote by £h the union of 
edges e (or faces for d> 3) of the mesh and the union of internal edges by £j^. Define P as 



(2.1) 
(2.2) 



eAu + b{x) ■ Vn -|- c{x)u = f{x) 

u = 
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on dCl 
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Figure 2.1: An example of a cdG decomposition. 



the union of boundary edges, i.e., those lying in d^l. The diameter of an element E €z Th 
is denoted He and h = maxE^Th ^E- We also denote the mesh function hy He, thus for 
X £ E we let hE^x) be equal to the diameter of E. We only consider meshes where h < 1. 
Define he := min(/i^- , hE+) for e G with e = E^ n for E'^, E^ G T^. The mesh is 
shape regular so there exists C > such that for all E we have < ~^^E-) — Ch^, 

and C > 1 such that C~^hE- < < CHe-- 

We define by Cl- decomposition the splitting of into two regions ^IcG and OdG such 
that for the closure Q = QqQ U J7dG) and we define by Th- decomposition the splitting of 
Th into two sub- meshes TqG and 7dG such that TcG C ilcC and Tdc '■= Th\ TcG- By 
abuse of language, we denoted here by %g not just the sub-mesh but also the region it 
occupies. Define Fcg (resp. T^g) to be the intersection of T with TcG (resp. Tdc)- Define 
J := TcG 1^ TdG and by convention we say that the edges lying in J are only part of 
the discontinuous Galerkin skeleton <fdG; the union of faces in TdC; and not part of the 
continuous Galerkin skeleton defined by ^cG '■= £h\ ^dG- 



In Figure 2.1 we illustrate a splitting for a problem where Q = (0, l)'^ and the solution 
exhibits layers at x = 1 and y = 1. The ^-decomposition is labelled, with the demarcation 
between the QcG and fidG regions given by a dashed line. A 7/i-decomposition is shown 
with the 7dG region shaded and the edges in J marked with a heavy line. 

Given a generic scalar field z/ : — )• M, that may be discontinuous across an edge 
e = E^ n E^ for E^, E^ G Th, we set := i^\e±, the interior trace on E^ and similarly 
define = tIe^ for a generic vector field r : $7 — t- M'^. Define the average and jump for 
a generic scalar as 

{M} ■=li^^ + ^~), M:=u+n+ + u-n-, on e e 

and for a generic vector field as 

{{t}} := 1{t+ + T-), H := r+ • n+ + T- • n-, on e € 

where is the unit outward pointing normal from E^ on e. For e G F the definitions 
become 

lu} := zv, M := un, {{r}} := r, on e G F. 
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We assume that the sign of 6 • n is the same for every x ^ e. Given a vector b denote 
the inflow and outflow boundaries of 0, by 

r"' ■.= {x£dn:b-n< 0}, r°"* ■.= {xedn:b-n>0} 

and for an element 

d'^'E :={x£dE:b-n<0}, d°'^^E := {x e dE -.b ■ n> Qi}. 

On an edge e, we denote by the trace of a function i/, taken from the element which 
contains e in its inflow boundary. 

Definition 2.3. The cG space is defined by 

V,G:={veH\n):yEeTh,v\EeF''} 

and the dG space by 

(2.4) V^a ■■= G : G Th,v\E S P'^} 

where is the space of polynomials of either total (or on quadrilateral meshes partial) 
degree at most k supported on E. 

We seek a finite element space V^^q that lies between the cG and dG spaces in the 
sense that 

(2.5) Kg C y,dG C T/dc- 

We construct this space by using continuous shape functions away from any boundary or 
internal layers present in the solution to (2.1 ) and discontinuous functions in the region of 
such layers. 

Definition 2.6. The cdG space is defined by 

(2.7) V,,a ■■= {v e L\n) : yE G v\e g p'; v\r,^a = 0; ^\%g ^ C^Tcg)}- 

Throughout we use the same polynomial degree k for V^q, V^q and V^dC- 
Let X be the characteristic function on Tdc, i-e., that defined by 



(2.8) X 
Then define 



I X £ TdG, 
xgTcg- 



^cdcCTdG) := {xv.ve Kdcl 

and 

KdG(TcG) := {(1 -X)v:ve V,^g}. 
In order to ensure a unique solution we make the following assumption. 

Assumption 2.9. We assume 

(2.10) r{x) := c{x) - ^V- b{x) >p>0 Vx £ n, 

for some p G M. 
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Definition 2.11. We define the local mesh Peclet number to be ||6||2,oo(£;)/i_B/(2e), see 

We consider meshes in the pre-asymptotic regime by making the following assumption. 

Assumption 2.12. We assume that for e = emax o-nd every E ^ Th the local mesh Peclet 
number is greater than \/hE. Moreover, we require max{||/i£;/6||^oo(-7-^^'), /i^;} < 1. 

As a consequence we have 
(2.13) e < emax < J min /i^/^||6||ioc(n). 

This assumption, for a fixed b, restricts the refinement of the triangulation for a given e. 
If we allowed /i — t- for fixed e > any layers would be resolved by the mesh and in the 
limit we would not see the non-physical oscillations associated with the cG approximation. 
We return to this question in Remark |6.6[ 

To characterise admissible 0-decompositions of the mesh we introduce the reduced 
problem: 

b ■ Vun + cnn = f on Q, 

(2.14) 

^ ' uo = onr°. 



Further, we define Ue := u — no, where u is the solution to the ADR problem (2.1). 

The ri-decomposition is chosen such that and uq have additional regularity on 
OcG- In general we do not expect that uq € i?^(0), even if we place higher regularity 
requirements on /, see, e.g., [3] and the references therein. 

Assumption 2.15. The set ^cG C 0, is chosen such that uq G H'^{^cg) o,nd \\us\\h'^{q^^q-) 
is hounded uniformly in e, that is for every < e < ffmax 

(2-16) \W\\m{n,G) ~ 1- 

2.1 Decoupled and Standard Formulations 

We define by V/^ the elementwise gradient operator. We discretize the advection term by 
(2.17) Ba{w,w):= Y [ {b-Vhw)wdx- Y [ b ■ Iwj w''' ds - ^ f{b-n)wwds, 

and the reaction term by 



(2.18) Br{w,w) := / cwwdx. 



E<^Tk 

The advection and reaction parts will frequently occur together and so for brevity we also 
define Bar{w, w) := Ba{w, w) + Br{w, w). 

For the diffusion term, besides the standard symmetric interior penalty method we 
also present a modified scheme, which we call the decoupled method. We refer to [1] for a 
comprehensive study of the interior penalty method and alternative discretizations. For 
w^w £ V^Q the decoupled method is defined by 



Bd{w,w) := Y / ^hW-Vhwdx 

+ E t 'yK^M-M-{iVHw]]-M + lVhw]]-M)ds, 



(2.19) ^^^^ 



ee£h\J 
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yielding 
(2.20) 



Be{w, w) := eBdiw, w) + Ba{w, w) + Briw, w). 



Here cr > is a discontinuity penalization parameter. With this formulation there is no 
control on the fluxes across J (hence the name decoupled). We recall that the standard 
symmetric Interior Penalty (IP) method is given by 



Bd{w, w) := Bd{w, w) + 

eeJ 



ah-'lwj ■ H - iiVhw}} ■ H + {{V,^^}} • H) ds, 



and 



Be{w, w) := eBdiw, w) + Ba{w, w) + Br{w, w). 



When restricted to the cdG space, the decoupled and standard IP forms become the 
bilinear form for the standard cG method on the continuous region. We first first analyse 
the stability properties of the decoupled formulation and infer stability for the standard 
IP method from a perturbation argument. 

We introduce the following mesh dependent norm for w £ y^G- 



(2.21) 
where 

and 



w 



|2 ._ 



H\d+\Mar^ 



Ml 



\W\ 



,V2...ll2 



w 



lL2{e)' 



with r defined in (2.10). 

Recall that for the symmetric interior penalty method the parameter a is selected 
independently of the Th such that Bd is positive definite with a coercivity constant which 
is also independent of Th-We assume that a is such that for all w G V^qI 



(2.22) 



\\{{VhW} • HIIli(^,) < h\'^hW\\L2in)\\V^eM\\LH£,)- 



Then, by Young's inequality, we have 

1„ 



w\\d < Bdiw,w) 



We adopt for Bd the same a as for Bd- 

We introduce a projection operator following the presentation of [2] which allows us 
to consider non-constant b. For polynomial degree fc > consider the L^-orthogonal 
projection IIj, : L^(r2) — )• F^^Q(7dG) defined by 

(2.23) / U^{v)wdx= / vwdx Vtt; e V^dcCTdc)- 
Jn Jn 

In particular 11-^(^)17;^, = 0. Furthermore, for all elements E of the mesh 

(2.24) \\Uj,iv)\\L2(^E) < \\v\\mE) yv E L\E). 
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As n-p(w) G KdG('^G) have for all E £ Th the inverse inequality 

(2.25) \U^{b ■ VhV)\HiiE) < hE^\\UT,ib ■ VhV)\\L2^E), 

and using a trace inequahty we have 

(2.26) ^ E h-E'm^{b-VHv)\\l,^Ey 



Define the streamline norm by 



(2.27) 

where te is defined by 
(2.28) 



bill := 111^111'+ rEmv{b-Vhv)\\l2 



Te '■= T mm 



hE hi 



\l'^(E) ^ 



and r is a positive number at our disposal. 



Definition 2.29. A decoupled cdG approximation to (2.1) is defined as Vh G ^cdG 
fying 



(2.30) 



Definition 2.31. A cdG approximation to (2.1) is defined as G V^^q satisfying 

(2.32) Be{vh,v)= [ fvdx \/vGV,^G- 

Jn 

We require that b points on J non-characteristically from TcG to TdG- 
Assumption 2.33. The Th- decomposition is such that for every e G J 



(2.34) 



^{bix)-n^)\e > emax^ 



Va; G e, 



where represents the unit normal pointing from TcG to TdG- 



Observe that the scaling between e and h mirrors that of Assumption 2.12 
Theorem 2.35. On V^q the bilinear forms Be and Be are coercive with respect to \\\w\\ 



(2.36) 



1, 



\vu\\\'^ < Beiw,w), 



1, 



IwlW"^ < Be{w,w), WGV^q. 



Proof. For the advection and reaction terms using integration by parts we have 



(2.37) Bariw, w) = ||rV2u;||i.(^) + E / 9 



For the diffusion term it follows from Equation 2.22 and Young's inequality that 
Bd{w,w) + y^^M-Mds> -||u;||^, Bd{w,w)>-\\w\\l 



Combing the last inequality with (2.37), the result now follows with Assumption 2.33 □ 



It follows that Vh and Vh exist and are unique. The following final assumption permits 
the use of an inverse inequality on the continuous Galerkin region. 



Assumption 2.38. The mesh TcG is quasi-uniform. 

It is convenient to denote the mesh-size on TqG by h-j-^^ 



\h 



£^lli°°(rcG)- 



The main result of this work is Theorem 6.2 , which states that the cdG approximation 



is stable in the streamline diffusion norm whenever Assumptions 2.9, 2.12, 2.15, 2.33 and 
12.381 are satisfied. 

In order to prove this result we establish first two separate stability bounds: We define 
Uej-uo £ KdG condition that for all v G V^cdG 



(2.39) 
(2.40) 



Observe that by linearity of the decoupled cdG method we have Ve + vq = Vh- In Section 
[3] we derive a bound for the decoupled cdG approximation Vg to on TcG] in Section |4] 
we obtain a bound for the decoupled cdG approximation vq to uq on TcG- In Section [5] 
we establish an inf-sup condition with streamline control on 7dG • Finally, in Section [6] we 
combine these results to show stability of the decoupled and then of the standard cdG 
approximation. 



3 Bounds on the Component on Tcq 

We introduce the projection operator of Scott and Zhang, [T^] |12l Section 1.6.2]. 

Lemma 3.1 (Scott-Zhang Projection). The Scott-Zhang operator SZ^ : W^'P{Q,) — t- V^q 
is a mapping with the following properties: For I > ^ there exists a Cgz > such that for 
all < m < min(l, /) 

(3.2) \\SZh{v)\\Hrn.^r.^)<CMmi%a) ^vGH'iTcG), 

and provided I < k -\- 1 for all E G T^c and < m < I we have the approximation 

(3.3) 11^; - SZh{v)\\H^iE) < Csz/i^^blHHAs) V?; G hHAe). 

where is the node patch of E, i.e., the set of cells in TcG sharing at least one vertex 
with E. 

Theorem 3.4. The decoupled cdG approximation is stable on the TcG region in the 
sense that 

(3.5) \\ve\\mi%G) ^ 1- 

Proof. We pick the auxiliary solution G KdG follows: On T^Gj define to be the 
Scott-Zhang projection of Us] and on 7dG to be the dG approximation with boundary 
conditions given by SZh{ue) on J and on TdC) i-e., 

VA = SZh{ue) on TcG, 
SeivA, v) = Be{Ue, v) \/v G KdGCTdc)- 
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Set 1] := Us — Vj[ and ^ := Vj\^ — Ve, so 1] + ^ 



■Ve- Notice that S V^^q. The Galerkin 



orthogonahty expressed by (2.39) and Theorem 2.35 give 

(3.6) \\m\\<m,i) = -Be{ri,0 = -Be{ri,i-X0, 



where x is defined in (2.8). Note that C ~ is continuous except on J where \^ — x^l 



■ rf and — = 5^''' '^i^^re the superscript C indicates the trace taken from the 
continuous Galerkin side of J. 

We examine each term of in turn. For the diffusion parts we use Young's inequahty 



For the advection term we use Assumption 2.33 which ensures that flux terms on J are 
zero as the upwind value of — vanishes. With Young's inequality we have 



where p is defined in (2.10). Finally for the reaction term 

4 ll„l|2 



g)- 



Using the previous three results, (3.6), the definition of the norm (|2.21|), and Lemma 3.1 
we gather ^ terms on the left hand side to show, with /it^^ 



4|l„l|2 |U||2 



(3.7) 



where in the final step we have used (2.16). As /? > we may use (3.7) and an inverse 
inequality to show 



(3.8) 

Assumption 



<l. 



2.15 



and ([3^ give ll^^lli,,.^ ^ < 1. 



□ 



4 Bounds on the vq Component on TcG 

We now pick the auxiliary solution as follows: On 7dG let it be uq and on 7dG be the 
dG approximation to uq with boundary conditions given by uq on F^g U J, i.e.. 



(4.1) 
(4.2) 



VA = uo on TcG, 



Lemma 4.3. We have for all v G V^^q that Be{vA,v) = Be{vo,v). 
Proof. Fix V £ F^^q. Then using ( |2.40[ ) 

Be{vo, v) = Beiuo,v) = Bs{uo, V - Xv) + BeiuQ, X^) 



where x is defined in (2.8). Observe that Bs{uo,xv) = BsivAjX"^) by (4.2). Notice that 
v—xv and uq are continuous on TcG- Recall that no integral over J appears in the definition 
of Bd- For Ba{vo, v — y f), th e integral over J vanishes since the value of {v — x^^)™ is zero 
because of Assumption 2.33 Therefore ;S£(?;o, f—x^) = i3£{uo,v—xv) = Be{vA,v—xv)- □ 
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Lemma 4.4. We have \\vo\\h'^{%g) 
Proof. Define Vj^ to be 



< 1. 



_ j SZh{uo) on TcG, 



on 



TdGi 



and let rj := v_a — v-j^, ^ := v-j^ — vq. With these definitions r] + ^ = v_a — vq, tjIj-^^ = and 
^ and r] are continuous on TcG- Then using Lemma 4.3 we have 

l\ur<m,c) = -Bdr^,o 



'ds. 



Due to Assumption 2.33 



- / eVhV ■ V/,e + {b ■ Vhv)C + C7]^dx+ [ b- MC 

have = the trace from the dG side of J, and [77]] = rj'^rf , 
' ' " We split each of the terms using Young's 



1 



ssumpiion z.oo we nave q, = i; , tne 
the trace and normal from the cG side of J 
inequality, giving 

_L ^11 l|2 II ||2 . P \ 

+ -nL--{n)\mm%G) + 1^1 

For the final term we note that ^ is a polynomial and so using 



(4.5) +^llc|li°°(f.)ll^lli^(reo) + ^ll^lli^(reG)+ / (^•^S')^''d.. 

i wc liULc uiiciL ^ ID a polynomial and so using Young's inequality and a 
trace and inverse inequality (with constant Cti) gives 



(4.6) 



j[b-n 7] )^ dS< — ||7? + — 



We combine (4.5 



and (4.6) to hide all terms of ^ under the norm on the left-hand side of 
(4.5). Using (3.3) for the terms of ij and a trace inequality gives 

Pmh^a) ^ infill' ^ (^^L +^r.o + ^L)ll^o|li.(r.o) ^ ^reGll^o|||.(r.o) 

and, by an inverse inequality, 11^11/^1(7-^) i$ 1- Now the result follows from the stability of 
the Scott-Zhang operator. □ 

5 Inf-Sup Condition 

The following theorem is an adaptation of related stability bounds in [5] and [2j to fit 
the above assumptions. Although the verification of the below inf-sup condition follows 
the overall structure in [5], we state it here in detail as the present analysis extends the 
scope to non-constant advection coefficients via the incorporation of II-p as [3]. Moreover, 
it deals with the modification of the bilinear form on J and it only has streamline control 
on the 7dG side. It is helpful to recall that H-pulT^Q = for any v. 

Theorem 5.1. There exists a positive constant A[s which is independent of h and £ but 



may depend on the polynomial degree, a, and the constants in (2.25) and (2.26) such that: 



(5.2) 



inf sup 



> A, 
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Proof. For an arbitrary v S V^^c define 

(5.3) v:=v + 'yvs, vs ■■= ^ TEllT^ib-Vhv), 



where 7 is a positive parameter at our disposal and te is defined in (2.28). Note that 
through the definition of 11^ we have v,vs G Theorem 5.1 is equivalent to showing 

the following two results: 

(5-4) \\v\\s < ll^^lls, 

(5.5) Be{v,v)>\\vfs. 

Consider first (|5.4|). We examine each term of ||fs'||c' in turn. We have 



(5.6) 



Yl ^\''s\h1{E) ^ Yl ^^E^\\^E'n^v{b-VhV)\\l2^E) 
EGTh EeTh 

< Y ^^EWUvib-'^hV^h^E) ^ \Ms- 



EeTh 



Also 



(5-7) \\r^^^vs\\l2(^Q) < IkllLoo(n) rl\\^vib-'^hv)\\l2(^E) ~ ll^lll- 

E(£Th 

For the terms on the edges we use ( |2.26[ ). This gives 

(5.8) ^|||^-nr/^MIli2(e) < MLo^in)r'EhE'mAb-Vhv)\\l2(^E) ^ Ms- 
Similarly, 

(5-9) E ^IIMIli^(e) ^ E 4^mv{b-V,v)\\l2^^^ < Ml 

The final term of the streamline norm gives 

Y rE\\^v{b-VhVs)\\h^E) < Y Mb-VhiTEUj,ib-Vhvml2^E) 

< Y 4\\b\\l^^E)hE''mvib-VHv)\\l2^E^ < 



Hi 



Combining the above results we have < 1 1 1 1 5- Using a triangle inequality we find 

\\v\\s < \\v\\s + 7\\vs\\s < CiT,a,j)\\v\\s, 



which concludes the proof of (5.4). 

To prove (5.5) first consider the advection and reaction terms of the norm. Using the 
linearity of Bar we have Bariv, v) = Bar{v, v) + ^Bari^, vs)- The second term is given by 

Bariv,Vs)=Y [ CviTEUj,{b-VhV)) + ib-VhV)iTEUj,{b-VhV))dx 

EeTh 

Y fb-lvjiTEUj,{b-Vhv)rds- Y [ib-n)v{TEUj,ib-Vhv))ds. 
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Using the properties of II-p given in (2.23) the second term above becomes 

^ / ib-Vhv)iTEll^ib-Vhv))dx= I TEllv{b-'^hv)Iijy{b-Vhv)dx 



(5.10) 



rEmv{b-'^hv)\\l2f^Ey 

E&Th 



Using Young's inequahty we have 



cv{TE'U-'D{b-yhv))dx <\\c\\L^{n) Y 2h\\l\E) + ifE\\^v{b-^hv)\\\2(^E) 



EeTh 



and, where C arises from a trace inequahty and the number of edges per element, 
Y ^b■lvjiTEUJ,ib■V^,v)rds- j{b-n)v{rEliv{b-^hv))ds 



< 



eer" 
2A ' 



EdTn 



In conclusion, using (|2.37), we have 

7l|c||L=o(f1) 



(5.11) 



El 

E&Th 



\l^(E) 



L2(e) 



■7 



E 

EdTh 



TE 



TET 

2A 



in. 



V 



{b-'^hv)\\l2^^y 



Recall that ||/i£;/b||L°°(7;jG) — ^ by Assumption 2.12, which imples te <^ for all G 7/i. 
For general all terms on the right-hand side of (5.11 ) are positive, provided A is large 
and 7 is small enou gh. 

Equation (2.22) ensures the continuity of Bd with respect to || • H^; thus 



(5.12) 



Bd{v,v) < Ci\\v\\d\\v\\d, 



for some Ci > 0. Recalling (5.6) and (5.9), it is clear that ||f5||d < ^211^11^, for some 
C2 > 0. Hence 

(5.13) Bd{v,v) = Bd{v,v) +'yBdiv,vs) > l\\v\\l - 7Ci\\v\\d\\vs\\d- 

Thus, if 7 < C1C2/8, then Bd{v,v) > which, combined with (5.11), gives (5.5). □ 



6 Stability of the Decoupled and Standard Approximations 

We saw that, under a set of suitable assumptions, the decoupled approximation satisfies 
the stability bounds: 



(6.1) 



\^h\\m{%G) ~ 1' 



\vh\\s < 



L2(Q). 



The first bound is a consequence of Theorem 3.4 and Lemma 4.4, the second of Theo- 
rem 5.1 So while one has streamline-diffusion stability on 7dG) even stronger bound 



is available on TcG under the aforementioned assumptions. We now derive a stability 
result for the cdG method. We require that the geometry of the interface J does not 
become significantly more complicated as the mesh is refined. More precisely, we require 
the boundedness of the trace operator. 
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Theorem 6.2. Suppose that the operator norm of the trace H^CTcg) ~^ ^'^{J) ^-^ hounded 
independently of h. Then, the cdG approximation Vh is stable in the sense that 

hTj\yhVh\\h(^T.^) + WvhWl < 1 + \\f\\h(ny 
Proof. Set Q := Vh — Vh- Using the coercivity of Be, Galerkin orthogonality and the norm 
of the trace H^{Tcg) — ^ L'^{J)i we have 



{V;.M-ICl + {v^C}-M 



< e[hr^jVhVh\\l2^n) ■ ^mWl^ij) + W^hCWh^n) 

he 



ea 



+ 



< 



ecr 



1/2 



e VcGl|Vh?^fe||i2(n) + , 1/2 1 
^ hJ 



Ph\\\L\J) 



+ ^IIIC1I 

Ale 



1/2 



L2(J) 



2.33 



for eahe 



and thus, using Assumption 
(6.3) \\C\f<ehrA^hVhfL^ 



+ 



ea 



-I7^IIMIIl2(j)- 



Dividing through by h^y^^ and using an inverse inequahty on /o||C||l2(£;) gives 
(6.4) ^r.Gl|V/.Clli2(7;3) < /irilllCr < ^l|V,^,||i2(^) + -f^||Mlli2(^). 

Using Assumptions 2.12 and [2.33 as well as (6.1), we bound each of the terms in (6.4). 
Using a triangle inequality on ||VC||^2(7^q) we conclude that 



< 1. 



To show that ||wh||s is bounded we establish an inf-sup condition for Be- Indeed, (5.4) 
may be used without change. It remains to transfer (5.5) to Be- The inequality (5.11) is 
still available as the discretesation of the lower-order terms did not change. Different is 
that we now use Bd{v,v) < Ci|||f||| • IH-OIH in place of (5.12), justified by Assumption 2.33 
Appealing to (5.6)-(5.9), one has \\\vs\\\ < C2\\\v 

(6.5) 



for some C2 > 0. Hence 
Bd{v,v) = Bd{v,v) +-fBdiv,vs) > l\\v\\l - tC'iIII^IIIIII^^sIII- 



For 7C1C2 small enough and A sufficiently large, 7 C'llHulH |||t;5| 
^Bar{v,v), using again the positivity of the terms in (5.11). 



is bounded by h\\v\ 



i + 

□ 



Remark 6.6. Due to Assumptions 2-12 and 2.33 the above stability hound is valid for 
the regime e < h^^'^\\b\\ioo(jiy For completeness we briefly outline here how Vh\\ ■ ||//i(n) 

stability of the cdG method is established if e > h^ • The stability proof is in this 

case easier because the mesh Peclet number is smaller. If is smooth or convex and the 
coefficients have sufficient regularity then u is in H'^{Q). Indeed for e > |/i^^||6||j;^cx)(q) 
the norm of u is uniformly bounded in e. Suppose that the mesh is quasi-uniform. By 
Ced's Lemma, with h := max£;/i£;, a standard argument gives 

\\u - VhWmin) ^ (1 + + h'^/'^)\\u\\H2(n), 

and thus Vh\\u - Vh\\H^n) ^ 1- 
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(a) Solution u given by ( 7.2 1 



(b) Solution given by ( 7.3 1 



Figure 7.1: Solution u and Ue for e = 10 



7 Numerical Experiment 

Let = (0, 1)^. We seek to solve 

(7.1) —eAu + {—X, —y) ■ Vu = —x — y 

with Dirichlet boundary conditions chosen such that the solution is given by 

Erf (x/^/2e) + Erf {y/V^s) 



(7.2) 



u{x, y) = X + y 



Erf (l/^/2e) 



where Erf is the error function defined by Erf(x) = Jo ^ "i^- ^'^^ < e ^ 1 this 
problem exhibits an exponential boundary layer along the outflow boundaries x = and 
y = of width 0{^/e). 

Away from the layers the boundary conditions on the inflow boundaries x = 1 and 
y = 1 are well approximated by 7/ — 1 and x — 1 respectively. The hyperbolic solution with 
these boundary conditions is given hy UQ{x,y) = x + y — 2. This gives 



(7.3) 



Ue{x,y) 



Erf {x/V2e) + Erf {y/V2e) 
Erf (l/^/2e) 



We plot (7.2) and (7.3) for e = 10 in Figure 7.1 note that away from the layers the 



solution Ug is close to zero. 

We let QcG be a set of the type (1 — 5, 1)^, 0<5<1. For each 6 the supremum 



sup 

eG{0,em£ 



\Ue\\H\n,a) 



is finite. The dependence of H^^elliTaj-j^^g) with respect to 5 and to e is illustrated in Figure 
E2l 



For this example c— |V-5 = 1, so Assumption 2.9 is satisfied. Further, we fix e = 10 
and consider uniform square meshes of edge length 2~^ so that Assumption 



6 



2.12 



^ ^ is also 

satisfied as the smallest local mesh Peclet number being 488.28. 

We define TcG = [1 — ^h, 1]^, where 6h = m2~^, m E {0, . . . , 32}. Note that, having 
fixed the mesh, 6h is a discrete parameter. The interface J is composed of the edges lying 
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10"'' 10"' 10"* 0.001 0.01 0. 



Figure 7.2: \\u^\\H'-i{n„G) different values of i5. 




Figure 7.3: Difference between cdG and dG solutions. 



on the lines y = 5^ for x > 5^ and x = 6^ for y > 5^. The smallest value of 6 • n is 5^ 
occuring on the edges containing the point {5h,Sh)- Thus, in this case, Assumption 2.33 
reads ^mh > and is satisfied for all m G {1,...,32} (note that this assumption 

is trivially satisfied when m = 0) for this choice of e and /i if o" < 200; in the shown 
computations a = 10. 

In Figure 7.3 we plot the L'^{Q) norm, ^/e weighted H^{Th) semi- norm, and norm 
of the jumps on (represented by [[•]]) for both the difference in the dG and cdG ap- 
proximations and the error in the cdG approximation. Note from Figure 7.3 that the 
difference in the approximations increases only very slowly until the final data points 
(where 7^g ~ Th)- When the continuous region covers the layer, non-physical oscillations 
pollute the approximation in function of Sh- 



in Table 7.1 we show the number of degrees of freedom (dofs) as the continuous region 



is increased. Reducing the degrees of freedom to approximately 30% of the dG method 
degrees of freedom results in only a very slight difference in the norm, thus showing that 
a considerable saving can be made without compromising stability. 

We finally remark that, at least for the example considered here, the choice of TcG 
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1 - 5 


dofs 


% of dG dofs 


y/£\\Vh{wh-Vh)\\H^{^n) 


dG 


4096 


100 


0.0 


8 X 2-^ 


3361 


82.1 


3.1157e-08 


16 X 


2417 


59.0 


6.8911e-08 


24 X 2-5 


2121 


51.8 


8.7544e-08 


30 X 2-5 


1457 


35.6 


1.7934e-07 


31 X 2-5 


1276 


31.2 


2.7896e-07 


cG 


1089 


26.6 


1.2444e-02 



Table 7.1: Degrees of freedom with e = 10 



leaving one layer of elements at the outflow boundary is optimal. Indeed, adding even a 
single element to the TcG region results in oscillations polluting the solution. For example, 
for 

TcG = [2-^ 1]' U ([0.5, 0.5 + 2-5] X [0, 2-5]), 
i.e., adding a sinlge element to TcG halfway along the x-axis, results in 



\\wh-vh\\L2i^^) = ^.7mQ X lo-^ 

\/i||V/,K - Vh)\\L\n) = 4.5008 x 10-^ 
a significant increase on the norms for TqG = 1]^- Notice that this choice of TcG 



violates Assumption 2.33 
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